Intro

Exploration of some different scenarios for financial burndown based on various assumptions around present value, annual expenditures, annual returns.

References

Inputs & Output

Inputs needed:

  • starting balance
  • annual drawdown
  • return on balance
  • number of years of drawdown period

Output is remaining balance at end of period.

Single Scenario

## loop through each year, new calc based on balance at end of each year to get year-by-year schedule
startbal <- -600000
draw <- 60000
return <- 0.04
yrs <- 10

start_bal <- startbal
sched_all_s <- data.frame()
for(y in 1:yrs){
  remain <- fv(r=return, n=1, pv=start_bal, pmt=draw, type=0)
  sched <- data.frame('year'=y,
                      'balance'=remain)
  sched_all_s <- bind_rows(sched_all_s, sched)
  start_bal <- remain*-1
}
  • Start: $600,000
  • Draw: $60,000 / yr
  • Return rate: 4%
  • Years: 10

At end of 10 years, you will have $167,780 left.

sched_all_s %>% ggplot(aes(x=as.factor(year), y=balance))+geom_line(group=1)+
  scale_y_continuous(labels=comma)+
  labs(title=paste0("$",format(remain, big.mark=",", digits=0), " left after ", yrs, " years."),
       y="", x="year from start")

Multiple Scenarios

## loop through each year, new calc based on balance at end of each year to get year-by-year schedule
scenarios <- tribble(
  ~startbal, ~draw, ~return, ~yrs,
  -600000, 60000, 0.04, 20,
  -700000, 60000, 0.03, 20,
  -700000, 60000, 0.02, 20,
  -800000, 60000, 0.02, 20,
  -800000, 60000, 0.03, 20,
  -1000000, 70000, 0.03, 20,
  -1000000, 60000, 0.04, 20,
  -800000, 60000, 0.04, 20
)
scenarios <- scenarios %>% mutate(
  scenario=paste0(as.character(startbal/1000*-1),"k-",as.character(draw/1000),"k-",as.character(return*100),"%-", as.character(yrs),"y")
)

sched_all_scenario <- data.frame()
for(s in 1:nrow(scenarios)){
  scen <- scenarios$scenario[s]
  start_bal <- scenarios$startbal[s]
  sched_all <- data.frame()
  for(y in 1:scenarios$yrs[s]){
    remain <- fv(r=scenarios$return[s], n=1, pv=start_bal, pmt=scenarios$draw[s], type=0)
    sched <- data.frame('scenario'=scen,
                        'year'=y,
                        'balance'=remain)
    sched_all <- bind_rows(sched_all, sched)
    start_bal <- remain*-1
  }
  sched_all_scenario <- bind_rows(sched_all_scenario, sched_all)
}
pbar <- sched_all_scenario %>% ggplot(aes(x=year, y=balance, fill=scenario))+geom_col(position = position_dodge())+
  scale_y_continuous(labels=comma)+
  labs(title=paste0("How much left after ", max(scenarios$yrs), " years in different scenarios"),
       y="", x="year from start")

ggplotly(pbar)
pline <- sched_all_scenario %>% ggplot(aes(x=year, y=balance, color=scenario))+geom_line()+
  geom_hline(yintercept=0)+
  scale_y_continuous(labels=comma)+
  labs(title=paste0("How much left after ", yrs, " years in different scenarios"),
       y="", x="year from start")

ggplotly(pline)

Simulations

Setup and Run Simulations

## number of simulations
nsims <- 100
## number of years to simulate
yrs <- 20
## parameters for starting, return, draw - exact numbers determined in loop
## min and max starting balance
start_min <- 700000
start_max <- 1000000
## return rate (real return, after inflation)
r_mean <- 0.03
r_sd <- 0.07
r_max <- 0.20
## draw 
draw_min <- 50000
draw_max <- 75000
## inflation rate - index draw against this
inflr <- 0.02

sim_all <- data.frame()

for(s in 1:nsims){
  sim <- s
  bal_start <- round(runif(n=1, min=start_min, max=start_max),0)
  bal_start_set <- bal_start
  bal_start_no_draw <- bal_start
  sim_sched <- data.frame()
  for(y in 1:yrs){
    rrate <- round(min(rnorm(n=1, mean=r_mean, sd=r_sd), r_max), 3)
    return <- ifelse(bal_start>0,round(bal_start*rrate, 0),0)
    draw_ni <- min(round(runif(n=1, min=draw_min, max=draw_max),0),bal_start+return)
    draw <- min(round(draw_ni*(1+inflr)^y,0),bal_start+return)
    bal_remain <- bal_start+return-draw
    bal_no_draw <- round(bal_start_no_draw*(1+rrate),0)
    sim_yr <- data.frame('sim'=sim,
                        'year'=y,
                        'start'=bal_start,
                        'rrate'=rrate,
                        'return'=return,
                        'draw_ni'=draw_ni,
                        'draw'=draw,
                        'balance'=bal_remain,
                        'bal_no_draw'=bal_no_draw)
    sim_sched <- bind_rows(sim_sched, sim_yr)
    bal_start <- bal_remain
    bal_start_no_draw <- bal_no_draw
  }
  sim_all<- bind_rows(sim_all, sim_sched)
}
sim_all$sim <- as.factor(sim_all$sim)
  • Number of simulations: 100
  • Number of years: 20
  • Starting balance: 700000 to 1000000 (uniform dist.)
  • Rate of return: average 0.03 ave, 0.07 std dev (normal dist. with max 0.2)
  • Annual draw: 50000 to 75000 (uniform dist.)

Visualize Progress Over Time

psim <- sim_all %>% ggplot(aes(x=year, y=balance, color=sim))+geom_line()+
  geom_hline(yintercept=0)+
  scale_y_continuous(labels=comma, expand=expansion(mult=c(0,0.01)))+
  theme(legend.position = 'none')+
  labs(title=paste0("How much left after ", yrs, " years in different simulations"),
       y="", x="year from start")

ggplotly(psim)

Summarize

Distribution of End Balances

chart_title <- "Distribution of End Balances at 20 yrs"
sim_all_20 <- sim_all %>% filter(year==20) 
sim_all_20 %>% ggplot(aes(x=balance))+geom_histogram()+
  geom_vline(xintercept=mean(sim_all_20$balance), linetype='dotted')+
  labs(title=chart_title)

What Really Counts: under water or not?

At 20 yrs
sim_all <- sim_all %>% mutate(
  status=ifelse(balance>0, "money","no_money")
)
sim_20 <- sim_all %>% filter(year==20) 
sim_20$status <- factor(sim_20$status, levels=c('no_money','money'))

#s_color <- c("light red", "dark green")
## get from RColorBrewer -> identify palette, select items; in this case, reverse order
s_color <- brewer.pal(n=3, name='Dark2')[2:1]

p1 <- sim_20 %>% ggplot(aes(x=status))+geom_bar()+
  scale_y_continuous(expand=expansion(add=c(0,5)))+
  labs(title="How often run out of money???", x="")

p2 <- sim_20 %>% ggplot(aes(x=1, fill=status))+geom_bar(position='fill')+
  scale_y_continuous(labels=percent, expand=expansion(add=c(0,0)))+
  theme(axis.text.x = element_blank())+
  scale_fill_manual(values=s_color)+
  #scale_fill_brewer(palette='Dark2')+
  labs(title="", x="", y="")

grid.arrange(p1, p2, nrow=1, widths=c(3,2))

Money vs No Money by Year
## summarize status results by year

sim_check_all <- data.frame()
for(y in 12:20){
  yr_check <- y
sim_check <- sim_all %>% filter(year==yr_check) %>% group_by(status) %>% summarize(status_count=n()) %>% mutate(yr=yr_check)
sim_check_all <- bind_rows(sim_check_all, sim_check)
}
## set factor levels for more control
sim_check_all$status <- factor(sim_check_all$status, levels=c('no_money','money')) 
## bar chart for split between money / no-money
p1 <- sim_check_all %>% ggplot(aes(x=as.factor(yr), y=status_count, fill=status))+geom_col(position='fill')+
  scale_y_continuous(labels=percent, expand=expansion(add=c(0,0)))+
  scale_fill_manual(values=s_color)+
  labs(x='years from start', y='')

## spread data to calc percentages of cases with money
sim_check_all_wide <- sim_check_all %>% pivot_wider(names_from=status, values_from=status_count) %>% mutate(
  status_pc=money/(money+no_money)
) 
## line chart showing percent with money each yr
p2 <- sim_check_all_wide %>% ggplot(aes(x=yr, y=status_pc))+geom_line()+
  scale_y_continuous(labels=percent, expand=c(0,0), limit=c(0,1))+
  scale_x_continuous(expand=c(0,0))+
  labs(x='years from start', y='% of simulations')

grid.arrange(p2, p1, nrow=1)

sim_check_all_wide %>% select(yr, status_pc)
## # A tibble: 9 x 2
##      yr status_pc
##   <int>     <dbl>
## 1    12      0.85
## 2    13      0.77
## 3    14      0.62
## 4    15      0.42
## 5    16      0.28
## 6    17      0.18
## 7    18      0.15
## 8    19      0.08
## 9    20      0.06
#sim_check_all_wide %>% select(yr, status_pc) %>% datatable()
sim_neg <- sim_all %>% filter(balance<=0)
sim_neg_yr <- sim_neg %>% group_by(sim) %>% summarize(yr=min(year))

sim_neg_yr %>% ggplot(aes(x=yr))+geom_bar()+
  labs(title="Distribution of Years when Balance = $0")

For those simulations that go below $0 with end balance.

chart_title <- "Dist. of Returns - Pos Sims"
hist1 <- sim_all %>% ggplot(aes(x=rrate))+geom_histogram()+
  geom_vline(xintercept=mean(sim_all$rrate), linetype='dotted')+
  ## geom_text preferred but geom_label used for better visibility
  geom_label(aes(label=mean(sim_all$rrate), x=mean(sim_all$rrate), y=+Inf), position='identity', vjust=1.2)+
  scale_y_continuous(expand=expansion(add=c(0,1)))+
  scale_x_continuous(labels=percent)+
  labs(title=chart_title)

mrate_all <- mean(sim_all$rrate)
sdrate_all <- sd(sim_all$rrate)

## calc returns when draw downs excluded
sim_rr <- sim_all %>% group_by(sim) %>% summarize(start_bal=first(start),
                                                  end_bal=last(bal_no_draw),
                                                  yrs=max(year)) %>%
  mutate(ttl_return=end_bal/start_bal-1,
         ttl_return_ave=ttl_return/yrs)

chart_title <- "Dist. Returns - Pos Sims (no draw)"
hist2 <- sim_rr %>% ggplot(aes(x=ttl_return_ave))+geom_histogram()+
  geom_vline(xintercept=mean(sim_rr$ttl_return_ave), linetype='dotted')+
  geom_label(aes(label=mean(sim_all$ttl_return_ave), x=mean(sim_all$ttl_return_ave), y=+Inf), position='identity', vjust=1.2)+
  scale_y_continuous(expand=expansion(add=c(0,1)))+
  scale_x_continuous(labels=percent)+
  labs(title=chart_title)

grid.arrange(hist1, hist2, nrow=1)

  • Mean rate of return (straight ave): 0.028617

  • Std dev rate of return (straight sd): 0.0698718

  • Mean total rate of return, no draws: 0.0373061

  • Std dev rate of return, no draws: 0.025411

chart_title <- "Returns over time"
sim_all %>% ggplot(aes(x=as.factor(year), y=return))+geom_boxplot()+
  geom_hline(yintercept = 0)+
  scale_y_continuous(labels=comma)+
  labs(title=chart_title)

chart_title <- "Draws over time"
sim_all %>% ggplot(aes(x=as.factor(year), y=draw))+geom_boxplot()+
  geom_hline(yintercept=mean(sim_all$draw), linetype='dashed')+
  scale_y_continuous(labels=comma)+
  labs(title=chart_title, x='year from start')

Evaluate Positive Outcomes

Identify Positive Outcomes

## identify sims where balance > 0 at year 20
sim_pos <- sim_all[sim_all$year==20 & sim_all$balance>0,]
sim_pos_rate <- nrow(sim_pos)/nsims
## filter full list for positive sims only
sim_pos_filter <- sim_pos %>% select(sim)
sim_pos_all <- left_join(sim_pos_filter, sim_all, by='sim')

## summarize individual sims
sim_pos_ind_smry <- sim_pos_all %>% group_by(sim) %>% summarize(
  start_bal=first(start),
  ave_rate=round(mean(rrate),3),
  ave_draw=round(mean(draw),0),
  end_bal=last(balance)
)

chart_title <- paste0("Ending Balance after ", yrs," yrs for positive sims")
sim_pos %>% ggplot(aes(x=reorder(sim, -balance), y=balance))+geom_col(position=position_dodge())+
  scale_y_continuous(labels=comma, expand=expansion(mult=c(0,0.05)))+
  labs(title=chart_title, x='sim')

Overall Sim Summary

sim_pos_smry <- sim_pos_ind_smry %>% summarize(
  ave_start=mean(start_bal),
  ave_rate=mean(ave_rate),
  ave_draw=mean(ave_draw),
  end_bal=mean(end_bal)
)
sim_pos_smry
## # A tibble: 1 x 4
##   ave_start ave_rate ave_draw end_bal
##       <dbl>    <dbl>    <dbl>   <dbl>
## 1   939863.   0.0563   76431.  328916

Individual Pos Sim Summary

sim_pos_ind_smry 
## # A tibble: 6 x 5
##   sim   start_bal ave_rate ave_draw end_bal
##   <fct>     <dbl>    <dbl>    <dbl>   <dbl>
## 1 17       989031    0.048    73978  272639
## 2 31       918633    0.069    75576  605805
## 3 38       974507    0.058    77199   45370
## 4 64       899359    0.061    76537  444082
## 5 69       992551    0.04     78739   59549
## 6 85       865098    0.062    76559  546051

Detailed Table

## DT to add table of all rows for each positive sim for inspection
datatable(sim_pos_all)

Positive Outcome Characteristics

What is distribution of starting balances within positive sims?
sim_pos_yr1 <- sim_pos_all %>% filter(year==1) 
sim_pos_yr1 %>% ggplot(aes(x=start))+geom_histogram()+
  scale_x_continuous(labels=comma)+
  scale_y_continuous(expand=c(0,0))+
  geom_vline(xintercept = mean(sim_pos_yr1$start), linetype='dotted', size=2)+
  geom_label(aes(label=prettyNum(round(mean(sim_pos_yr1$start),0), big.mark = ","), x=mean(sim_pos_yr1$start),
                ## y set to +Inf which goes to top
            y=+Inf,
            hjust=1.2,
            vjust=1.5))

What is distribution in returns within sims?

chart_title <- 'Distributions in Annual Returns across positive sims'
sim_pos_all %>% ggplot(aes(x=sim, y=rrate))+geom_boxplot()+
  geom_hline(yintercept = 0)+
  geom_hline(yintercept = mean(sim_pos_all$rrate), linetype='dotted', size=1.5)+
  geom_text(aes(label=paste0(round(mean(sim_pos_all$rrate),3)*100,"%"), x=0, y=mean(sim_pos_all$rrate)),
            hjust=0, vjust=-1)+
  scale_y_continuous(labels=percent)+
  labs(title=chart_title)

What is the distribution of draws?

Not inflation adjusted, for comparison with current dollar value.

sim_pos_all %>% ggplot(aes(x=draw_ni))+geom_histogram()+
  scale_x_continuous(labels=comma)+
  geom_vline(xintercept = mean(sim_pos_all$draw_ni), linetype='dotted', size=2)+
  geom_text(aes(label=prettyNum(round(mean(sim_pos_all$draw_ni),0), big.mark = ","), x=mean(sim_pos_all$draw_ni),
                ## y set to +Inf, which goes to top; requires vjust
                y=+Inf, hjust=1.2, vjust=1.5))

Is there correlation between starting balance & returns?

Do sims in the positive set have both high starting AND high returns?

sim_pos_ind_smry %>% ggplot(aes(x=start_bal, y=ave_rate))+geom_point()